Understanding the Response of Poly(ethylene glycol) diacrylate (PEGDA) Hydrogel Networks: A Statistical Mechanics-Based Framework

Thanks to many promising properties, including biocompatibility and the ability to experience large deformations, poly(ethylene glycol) diacrylate (PEGDA) hydrogels are excellent candidate materials for a wide range of applications. Interestingly, the polymerization of PEGDA leads to a network microstructure that is fundamentally different from that of the “classic” polymeric gels. Specifically, PEGDA hydrogels comprise PEG chains that are interconnected by multifunctional densely grafted rod-like polyacrylates (PAs), which serve as cross-linkers. In this work, we derive a microstructurally motivated model that captures the essential features which enable deformation in PEGDA hydrogels: (1) entropic elasticity of PEG chains, (2) deformation of PA rods, and (3) PA–PA interactions. Expressions for the energy-density functions and the stress associated with each of the three contributions are derived. The model demonstrates the microstructural evolution of the network during loading and reveals the role of key microscopic quantities. To validate the model, we fabricate and compress PEGDA hydrogel discs. The model is in excellent agreement with our experimental findings for a broad range of PEGDA compositions. Interestingly, we show that the response of PEGDA hydrogels with short PEG chains and long PA rods is governed by PA–PA interactions, whereas networks with longer PEG chains are dominated by entropy. To enable design, we employ the model to investigate the influence of key microstructural quantities, such as the length of the PEG and the PA chains, on the macroscopic properties and response. The findings from this work pave the way to the efficient design of PEGDA hydrogels with tunable properties and behaviors, which will enable the optimization of their performance in various applications.


INTRODUCTION
Hydrogels are soft materials with tunable mechanical and physical properties, which are governed by their chemical composition.The wide range of properties and versatility of hydrogels enables high innovation potential in material technology.The response of hydrogels under a mechanical loading is of particular interest for a wide range of applications, including tissue engineering, 1,2 drug delivery, 3−5 traction force microscopy, 6,7 and microfluidics. 8,9In these systems, poly-(ethylene glycol) (PEG) based hydrogels were used thanks to their stability under various conditions and biocompatibility.
−9 The step-growth approach, e.g. by crosslinking bifunctional PEG with multifunctional cross-linkers of functionality f ≥ 3, results in relatively uniform PEG gel networks with point-like junctions.Figure 2a illustrates such an "ideal" hydrogel network formed using point-like tetrafunctional junctions which are cross-linked by random-coil PEG chains.PEGDA photopolymerization, on the other hand, offers numerous advantages, including facile hydrogel formation as well as the ability to engineer complex multidimensional hydrogel structures. 12,13−16 Figure 2b illustrates a hydrogel of randomly oriented bottlebrush-like polymers with rod-like PA cores that are intertethered and connect randomcoil PEG chains (Figure 2c). 17,18ritical to interpreting the structural and mechanical properties of PEGDA networks is understanding that each PEGDA macromonomer is capable of acting as both a monomer that enables the extension of PEG-based chains, as well as a cross-linking agent that connects two separate acrylate (A) cores.Notably, both acrylate groups on a single PEGDA macromonomer polymerize independently.Therefore, as the free radical photopolymerization of PEGDA proceeds, acrylate groups added to a radical PA chain end may originate from unreacted PEGDA macromonomers, or from PEGDA macromonomers in which one of the two acrylate end groups has already reacted and is fixed to another PA chain.This process results in a highly interconnected network topology.Moreover, adding acrylate to the growing PA rods necessarily involves the addition of PEG chains, which are grafted to the PA core.−21 Accordingly, the propagation reactions throughout the system result in a three-dimensional network structure consisting of rod-like PA chains, which are randomly oriented in space and serve as cross-links that interconnect entropic PEG chains.In this context, we point out that PEGDA-based hydrogels resemble a swollen molecularly reinforced composite (cf.refs 22, 23) with stiff PA chains acting as reinforcer (similar to a fiber in a fiber-reinforced composite) in the swollen PEG matrix.The matrix-reinforcer coupling is perfect since each acrylate repeat unit is chemically linked to the PEG matrix.Topological defects such as loops and dangling ends carrying an unreacted acrylate group may be present in the network.However, it seems reasonable to assume such structural irregularities are negligible in view of the overall bottlebrush network structure that determines the mechanical properties and response of the highly cross-linked PEGDA hydrogel network.
From the literature data discussed above, it is clear that the structure of PEGDA-based hydrogels obtained by free radical polymerization is characterized by densely cross-linked bottlebrush polymers.In particular, this is concluded from SAXS and SANS data of a series of hydrogels made from PEGDA of different molecular weights, 18,24 and estimates of the molecular weight of the polyacrylate chain, the key structural component of bottlebrush polymers (backbone), by selective ester cleavage of PEGDA-based networks. 25o the best of our knowledge, such complex bottlebrush networks have not been explored through theoretical modeling of PEGDA.Most available modeling approaches employ theories which either assume point-like junctions of functionality f = 3 to f = 4 or assume a star-like structure for f > 4. 17,24,26−31 In the latter case, a structural model derived from experimental small angle neutron scattering data provided a systematic explanation of the linear elasticity and swelling of PEGDA hydrogels in which the PA backbone was assumed to be short relative to that of the PEG chains. 17However, this approach neglected considerable evidence that long bottlebrush structures are the key structural feature. 18,24,254][25][26]32 The fundamental difference between the classical networks and the network formation and topology of PEGDA has farreaching consequences on the mechanical properties and constitutive behavior of the hydrogel. Theim of this contribution is to better understand the structure−property relationships and the overall response of PEGDA hydrogels under mechanical loading.To this end, we derive a microstructurally motivated and energy-based model that captures the three essential features that enable deformation in PEGDA hydrogels: (1) the entropic contribution of the PEG chains, (2) the deformation of the PA rods, and (3) the mechanical PA−PA interactions.We demonstrate that in the limit of short PEG chains and long PA rods, the behavior is governed by PA−PA interactions, which we attribute to the mechanical restrictions due the high density grafting of PEG chains to the PA cores.In the case of long PEG chains and short PA rods, the hydrogel behaves entropically due to the deformation of the PEG chains.We also develop simple closed-form approximations for the expressions for the Young's modulus and the stress.The model captures the key features of mechanical response and the microstructural evolution of the network during loading, and we show that it agrees with experimental findings for a broad range of PEGDA compositions.

MODELING THE RESPONSE OF PEGDA NETWORKS
In developing our microstructurally motivated and energybased model, we consider a PEGDA hydrogel network comprising N PEG entropic PEG chains and N PA stiff densely grafted PA chains per unit undeformed swollen volume.Each PA chain comprises n PA acrylate repeat units and is connected to n PA PEG chains.Accordingly, we can estimate N PA ≈ 2N PEG /n PA .
To model the overall response of the PEGDA hydrogel, we denote the material points in a referential traction-free unswollen configuration by X.The network is placed in an aqueous bath and swells.In response to an external loading, the network deforms.The material points in the deformed state are denoted by x and we define the deformation gradient from the reference to the deformed state F = ∂x/∂X.The deformation can be written as F = F s F m , where F s = J 1/3 I, where I is the identity matrix, is associated with the swelling of the network and F m denotes the mechanical deformation due to the external loading.It is assumed that the mechanical loading is isochoric and accordingly det F m = 1.
The deformation of PEGDA networks is enabled by the entropic extension of the softer PEG chains, which transfer forces to the stiff PA chains, as illustrated in Figure 3.We assume that the total energy-density of the deformed gel can be viewed as the sum of three contributions: where ψ PEG is the entropic energy-density due to the deformation of the PEG chains, ψ PA is the elastic energydensity resulting from the bending of the polyacrylate chains, and ψ int is the energy-density due to the mechanical interactions between the PA chains.Accordingly, the true stress is pI where the σ PEG = 1/J(∂ψ PEG /∂F)F T , σ PA = 1/J(∂ψ PA /∂F)F T , and σ int = 1/J(∂ψ int /∂F)F T denote the stress associated with the PEG chains, the PA backbone, and the mechanical PA−PA interactions, respectively, and p is a work-less pressure-like term that stems from the incompressibility of the network and is determined from the boundary conditions.It is also convenient to define the first Piola-Kirchhoff stress per unit swollen undeformed area: In the following, expressions for the three energy contributions described in eq 1 are derived.For convenience, Table 1 summarizes the model parameters that will be used.
2.1.Energy-Density of PEG Chains.The PEG chains are modeled as freely jointed chains (FJCs) with n freely jointed repeat units of length l.In the reference dry state, the end-to- , where R n l = and R i ( ) are the end-to-end distance and the chain direction, respectively.−40 Following common practice, we consider a network in which all chains assume the average referential end-to-end distance.By employing the affine deformations assumption, the deformed end-to-end vector of the i-th chain is where r r r are the deformed endto-end distance and chain direction, respectively.We highlight that the deformed end-to-end distance depends on the orientation of the chain with respect to the loading direction.
The entropy of the i-th polymer chain is given by 35,41 S n const k n ( , ) ( ) ln sinh ( ) where k b is the Boltzmann constant, ρ (i) = r (i) /nl is the ratio between the end-to-end distance and the contour length of the PEG chain, and β is determined from the Langevin function ρ = coth β − 1/β.A useful approximation that will be used in this work is β = ρ(3 − ρ 2 )/(1 − ρ 2 ). 42The free energy of the ith chain is ψ PEG (i) = −TS c (i) , where T is the temperature, and accordingly the stress associated with the chain is The overall energy of the PEG chains in the network is determined via where the summation is carried over all chains.Henceforth, the notation ⟨•⟩ = ∑ i = 1 N • (i) /N denotes the average of the quantity •.The stress that develops in the PEG chains is where the stress associated with a chain σ c is given in eq 6.

Energy-Density of PA Rods.
As described in the introduction, the acrylate (A) groups in the PEGDA chains polymerize to form a PEG network with multifunctional cross-

Macromolecules
links in the form of densely grafted PA rods.Once the PEG chains stretch, the entropic forces work toward bending the PA rods, as illustrated in Figure 3.In the following, we model the elastic energy-density of the PA rods and their contribution to the overall macroscopic response.To this end, we begin by examining a single PA rod that is aligned along a given direction and derive an expression for its energy.Next, we sum over the energies of all PA rods along that direction.To determine the overall energy, an additional integration process over all PA direction is carried out.Following the chemistry of the PEGDA network, the PA chains are treated as stiff rods that undergo bending due to the entropic forces exerted by the PEG chains, which are grafted to them. the plane perpendicular to the PA axis (see Figure 4).The average PA rod connects n PA PEG chains that are equally spaced along the PA with a distance l A = 0.25 nm (corresponding to the length of a repeat unit A) such that the overall length is L PA = n PA l A . 15,43Next, we assume that (1) the PEG chains are perpendicularly connected to the PA such the entropic forces from the PEG chains are uniformly distributed along the PA chain, and (3) the effect of dangling ends is negligible.
Following assumption (1), the end-to-end vectors of the PEG chains along the PA backbone are found on the plane perpendicular to N

PA ( )
. To characterize this, consider one of the PAs aligned along N

PA ( )
. We mark this specific PA rod by an index K, where 1 ≤ K ≤ N PA (α) .This rod connects n PA PEG chains, and we denote the end-to-end direction of the k-th where 0 ≤ φ < 2π is the angle between the end-to-end vector R K k ( , , ) and Y PA ( ) .Note that α, K, and 1 ≤ k ≤ n PA denote the direction of the PA rod, the specific configuration of the interconnected PEG chains along the PA rod, and the location of the PEG chain along the PA rod, respectively (see Figure 4).Upon the application of an external force, the PEG chains deform affinely (see eq 4) and rotate the PA rods.We denote the direction of the PA in the deformed configuration by n .
PA ( ) Assuming that the orientation of the PEG chains r̂with respect to the PA rods remains such that r n PA ( ) , we find that The entropic forces from the PEG chains result in the bending of the PA chain and are schematically illustrated in Figure 4.The force applied by the k-th PEG chain on the K-th Accordingly, the moment at a distance γl A ≤ ξ<(γ+1)l A , where 0 ≤ γ ≤ n PA is an integer, can be written as ) The full derivation is given in section S1 of the Supporting Information.Next, the overall energy of the K-th PA rod can be computed via Here, l p = EI/k b T is the persistence length of a densely grafted PA rod, with EI denoting the bending stiffness.In addition, note that since the integration is along the length of the rod, the energy density ψ r (α,K) does not contain the index k.To determine the energy of all PA backbones aligned along the direction , we follow the assumption that the PEG chains are uniformly distributed and integrate The integral in the last equality denotes the average energy of a chain aligned along the direction N PA ( ) and the resulting energy ψ r (α) does not contain the index K.
The average stress associated with the PA rods along the N PA ( ) direction can be determined from eq 14 via where and A = (2βdβ/dr)/(ρn), with dβ/dr = β 2 sinh 2 β/(sinh 2 β − β 2 ).Lastly, the energy-density of all PA chains is given by summing over all directions N PA (α) : Accordingly, the contribution of the PA rods to the overall stress is We point out that the stiffness of the PA rods stems from the collective effect of the densely grafted side-chains.Accordingly, two remarks are underscored.First, as opposed to the classical Macromolecules theories, the unique composition of PEGDA hydrogels results in a dependence between the entropic PEG chains and the PA rods.Specifically, adding monomers to the PA rod requires the addition of PEG chains.Accordingly, and as discussed in previous works, 44,45 the classical definition of persistence length to bottlebrush polymer backbone is somewhat problematic and the contributions of the side-chains (or PEG in this case) must be taken into account.Due to steric hindrance, excluded volume, and intermolecular interactions, the addition of A units (i.e., an increase in n PA ) may change the persistence length l p of a PA rod.−48 Second, in networks with short PA rods that interconnect only a small number of PEG chains, the bending energy stored in the PA is negligible in comparison with the entropic energy of the PEG chains.Therefore, in such networks the contribution of the PA rods can be neglected.

Interaction Energy-Density of PA Backbone.
As opposed to the classical networks, during the formation of PEGDA hydrogels, A groups connect to form a network of PA rods that interconnect to one another through PEG chains.As a result, the PA rods can interact mechanically, thereby significantly changing the underlying mechanics of the system.In this context, we emphasize that the PA rods are not likely to physically come into contact.Rather, mechanical interactions include steric restrictions of the bottlebrush polymers due to their shape anisotropy.These restrictions are governed by two main factors: (1) the rotational freedom of the PA rods due to their hydrodynamic volume and (2) the ratio between the lengths of the PA rods and the PEG chains.The former is caused by the entropic forces from the PEG chains which act on the PA rods.To account for the latter, we examine the ratio L PEG /L PA to account for these interactions (see Figure 5).
Recall that L PEG = nl and L PA = n PA l A are the contour lengths of the PEG chain and the PA backbone, respectively.In the case of L PEG /L PA ≪ 1 the PEG chains are significantly shorter than the PA backbone.Due to the stochastic nature of network formation, a PA rod can be interconnected to many other PA rods through mediating PEG chains, thereby increasing the density of PA rods per unit volume.As a result, mechanical interactions between PA rods is highly probable and is expected to affect the overall response.Since the PA rods are stiffer than the PEG chains, such interactions result in the stiffening of the network.If L PEG /L PA ≈ 1, the likelihood of mechanical interactions is lower, but still exists and results in a slight stiffening.In networks with long PEG chains such that L PEG /L PA ≫ 1, the PA backbones are far from each other and any PA−PA interactions can be neglected.
To model this effect, we propose the simple first-order interaction energy where I 1 =Tr(F T F) is the first invariant and α denotes the maximum stiffness contribution of the interactions in the limit L PEG /L PA → 0 (or in the case in which the PA rods are very close to each other).It is emphasized that in the limit L PEG /L PA → ∞, the PA backbones do not interact and the energetic contribution ψ int → 0. The interactions between the PA backbones depends on their density, and therefore it is reasonable to assume that α ∼ N PA .
The stress due to these intermolecular interactions is determined from eq 19 via

Integration from the Local Chains to the Network Level.
To integrate from the chain to the network level, we employ the microsphere technique.This numerical method enables one to determine an average stress contribution by integrating over the stress associated with the individual chains.−54 Consider a network with chains that are initially uniformly distributed and randomly oriented.The end-to-end vector directions R i ( ) can be represented as unit vectors that begin at the center and end at the surface of a unit sphere.The directional averaging of a quantity • over the surface of the unit sphere can be approximated via the discrete summation where the index i = 1, 2, ..., m refers to the value of • (i) along the i-th direction and w (i) are non-negative weight functions constrained by ∑ i = 1 m w (i) = 1.It is important to note that due to the uniform distribution and random orientations or the PEG chains, the averages R 0 = and R R I 1/3 = .In this work, we follow the findings of ref 49, which demonstrated that a specific choice of 42 integration directions provides sufficient accuracy in the computation of averages.The integration directions and the corresponding weights are given in Table 1 of ref 49.2.5.Closed-Form Approximation for the Stress in PEGDA Hydrogels.The full model involves numerical integrations and the use of the microsphere approach, which leads to long computational times and makes it harder to explicitly determine pertinent quantities such as stiffness and overall response.To gain deeper insight into the overall microstructural evolution, simpler and efficient closed-form solutions are required.In the following, we derive a closedform approximation for the stress that develops in PEGDA hydrogels.To this end, we consider a network with PEG chains and PA rods that are initially randomly oriented and uniformly distributed.
The stress that develops in the PEG chains can be approximated by considering the ratio ρ = r/nl ≪ 1. Accordingly, the PEG chains can be treated as Gaussian chains and the Langevin function can be approximated via a first order Taylor series as β ≈ 3ρ.The resulting stress associated with the deformation of the PEG chains is Next, we approximate the stress that develops in the PA rods.To this end, we take a second order Taylor series in ρ around zero in eq 15, substitute the resulting expression into eq 18, and compute the averages (see section S2 in the Supporting Information).This approximation yields and I 1 = Tr(F m F m T ) is the first invariant.The interaction stress given in eq 20 is a first order approximation, and accordingly the total stress is estimated by substituting eqs 22, 23, and 20 into eq 2. The first Piola-Kirchhoff stress is determined from eq 3.

COMPARISON TO EXPERIMENTS
3.1.Experimental Protocol.Hydrogel disks with high geometric and compositional uniformity were fabricated under oxygen-free argon atmosphere by using a 3 mm thick PDMS plate with cylindrical bores of 6 mm in diameter sealed by a light-transparent coverslip at both sides.After filling the bottom-sealed cylindrical cavities with deionized aqueous solution of PEGDA containing lithium phenyl (2,4,6-trimethylbenzoyl) phosphinate (LAP) as photoinitiator (4.95 mM) and sealing the upper opening, gel samples were obtained by free radical photopolymerization (5 min irradiation at 370 nm with 4 mW/cm 2 from both sides).In the gel syntheses, PEGDA concentrations c in the preparation state were above the critical overlap concentration c* of the corresponding PEG chain (recall that c* is the inverse of the intrinsic viscosity of PEG 55 ).By setting c > c*, we avoid significant inhomogeneities in the PEGDA hydrogel network. 24,56,57Depending on the numberaveraged molecular weight M n of PEGDA, determined by gel permeation chromatograhpy (GPC with PEG calibration) to be M n = 1.07 kDa, M n = 2.03 kDa, and M n = 5.85 kDa, the concentration c varied between ∼8−26%wt.For of given M n , three concentrations c were used.
To determine the volumetric deformation J (or the water content in the equilibrium swollen state of the PEGDA hydrogel discs), we weighed the hydrogel in the fully swollen and the dehydrated states.The volumetric deformation J, defined as the ratio between the swollen and the dry volumes of the gel, can be computed via where ρ p ≈ 1.18 g/cm 3 and ρ w ≈ 1 g/cm 3 are the densities of PEGDA and water, respectively, and W s and W d are the weights of the fully swollen gel and the dry polymer, respectively.We point out that different works report densities in the range of ∼1.07−1.18g/cm 3 of PEGDA. 26,27,29,58niaxial compression tests were carried out on the fully swollen samples using a texture analyzer (TA.XT Plus Connect, Stable Micro Systems) with a 50 mm-diameter metal plate fixture attached to a vertically movable load cell.Before each compression test, to prevent the hydrogel from sticking to fixtures, a silicone oil coating was applied to the top and the bottom surfaces of the hydrogel as well as the metal fixtures of the testing device.Compression tests were executed at a rate of 0.01 mm/s up to rupture. Figure 6 shows images from the beginning, the middle, and toward the end of a uniaxial compression test on a cylindrical sample with molecular weight M n = 2.03 kDa and concentration c = 18.4%wt/wt.A video of the experiment is given in movie SV1 with a description in section S3 of the Supporting Information.To ensure reproducibility, 3 specimens of each molecular weight and concentration were tested.The results are highly reproducible.

Model Fit. 3.2.1. Approximation -Uniaxial Tension/ Compression.
In the following we define a global deformation gradient x y z , , { } and examine the case of uniaxial compression.Accordingly, the deformation gradient due to a mechanical loading can be written as To determine the Piola-Kirchhoff stress, we employ the boundary conditions y Py z Pz 0 • = • = to obtain the pressure term p.In the limit of small to moderate deformation, we employ the approximations in eqs 22 and 23 to obtain where is the Young's modulus of the PEGDA network.Here, the quantity η is given in eq 24.In the case of L PA /L PEG ≪ 1, it can be shown that the second and the third terms in eq 28 are negligible in comparison to the first and accordingly the classical expression for the Young's modulus E ≈ 3N PEG k b T/ J 1/3 is recovered. 34,35.2.2.Model Fit.To validate the model, we compare its predictions to experimental findings.To this end, we consider three molecular weights M n of PEGDA and three concentrations c for each molecular weight.The temperature is T = 300 K and we assume an interaction stiffness α = N PA .In addition, the acrylate length l A , the persistence length l p of the densely grafted PA rod, and the length of a freely jointed repeat unit in the PEG chain l are approximated according to existing experimental data.Specifically, l A = 0.25 nm, 15,43 l p = 17.5 nm, 20,59,60 and l = 0.3 nm. 61The number of repeat units in the average PEG chain is estimated via n ∼ M n /M EG , where M n and M EG ≈ 44g/mol are the molecular weights of the PEG in PEGDA (i.e., PEGDA without the end groups (A)), and an EG monomer repeat unit, respectively.The contour length of a PA rod L PA and the chain-density of PEG N PEG were determined based on a fit to the experimental data.We comment that our experiments show that the dry PEGDA network is brittle, suggesting that the response is not entropic.For convenience, the parameters are listed in Table 2.
The fit of the model and the approximation in eq 27 to the experimental findings for PEGDA hydrogel networks with the molecular weights M n = 1.07 kDa, M n = 2.03 kDa, and M n = 5.85 kDa are shown in Figure 7a−c, respectively.The different

=
: ratio between initial end-to-end distance and contour length of PEG chain.concentrations are denoted in the plots.Here, the circle marks denote the experimental findings, the continuous curves the full model, and the dashed curve the closed-form approximation in eq 27.We find that the approximation provides a great fit for small to moderate stretches, with a maximum error of <17%.We find that as we increase the length of the PEG chains (i.e., increase n), the predicted length of the PA backbone decreases, resulting in a lower n PA .We believe this to be a consequence of steric hindrance and intermolecular interactions between the PEG chains during the gel synthesis.In addition, the model predicts that increasing the concentration of PEGDA while maintaining a constant molecular weight M n results in longer PA rods with the same PEG lengths (i.e., same n with larger n PA ).Our experimental measurements show that this is also accompanied by a reduction in the volumetric deformation J.Such trends are to be expected from the kinetics of PEGDA hydrogel network formation. 23,62nother interesting finding is that the term k b TN PEG increases with concentration c for a given molecular weight M n .However, comparison of this value under different molecular weights does not reveal a specific trend.These fits stem from the choice of α = N PA , which is set for all molecular weight (M n ) values.Specifically, in the case of high molecular weights, such as M n = 5.85, the overall behavior is governed by the entropic contribution of the PEG chains.In low molecular weights (i.e., M n 1.07 kDa), the stiffness may stem from two main sources: the bending energy and the mechanical interactions between the PA rods.The former is due to the high value of η, which results from the length of the rods.Given the proposed interaction energy and interaction stress (eqs 19 and 20) and the assumption that α = N PA , which denotes the maximum stiffness contribution of the mechanical PA−PA interactions in the limit L PEG /L PA → 0, we find that the interactions between the PA rods are negligible in the examined PEGDA hydrogels.Such interactions are expected to be significantly more dominant as L PEG /L PA → 0, and further investigations into the influence of the molecular weight M n and the network composition on the parameter α and the interaction stress is required.
To emphasize the influence of the different stress components, we plot the contribution of the stress components P PEG , P PA , and P int due to the deformation of the PEG network, the bending of the PA rods, and the PA−PA interactions, respectively, as a function of the stretch λ for hydrogel networks in Section S4, Figure S3 of the Supporting Information for three representative PEGDA networks.We show that in the case of M n = 2.03 kDa and c = 15.8%wt/wt, the PEG chains are shorter than the PA rods and accordingly the stress component P PA associated with the bending of the PA rods governs the response of the network.In networks characterized by M n = 5.85 kDa and the concentrations c = 8.0%wt/wt and c = 15.9%wt/wt the chains are significantly longer than the PA rods.However, the water content in the two PEGDA networks is significantly different, resulting in variations in the extensions of the chains.Chains that experience larger stretches significantly stiffen the overall network. 35,37In all cases the interaction energy is negligible, and is expected to be more dominant as the ratio L PEG /L PA → 0. However, we emphasize that a more rigorous experimental investigation into the interaction energy is required.
To further illustrate the influence of the bending of the stiff PA rods, we compare the model to the classical solutions from rubber elasticity of gels which consider point-like junctions 33,36,37 in Section S5 of the Supporting Information.It is shown that while all models capture the right trend, neglecting the contribution of the PA rods results in a softer network.Therefore, it is essential to account for the contribution of the PA rods in the analysis of PEGDA networks.
A comparison between the Young's moduli obtained from the experimental data and the proposed model is depicted in Figure 8.Here, the triangle, the square, and the circle marks denote the experimental data and the model prediction (eq 28) for M n = 1.07 kDa, M n = 2.03 kDa, and M n = 5.85 kDa, respectively, and the dashed line has a slope of 1 to enable comparison.Specifically, the distance of the points from the dashed line denotes the error in the model predictions.We find that the model is capable of capturing the stiffness with an error of <10% across a broad range of PEGDA molecular weights.It is emphasized that the stiffness in high molecular weights stems from the entropic contribution of the PEG chains, and therefore the first term in eq 28 is dominant.In lower molecular weights, the PEG chains are short and the PA rods are long, and therefore the stiffness is determined by the PA rods.Specifically, the second term in eq 28 is dominant.As previously described, the interaction term is negligible in the chosen hydrogels.

PROGRAMMING THE RESPONSE OF PEGDA NETWORKS
One important advantage of our microscopically motivated model is the ability to investigate the influence of microstructure on the mechanics of the hydrogel network.With the aim of tuning the properties and the response of PEGDA networks, in this section we explore the influence of different microstructural quantities.As above, we set T = 300 K, α = N PA , l A = 0.25 nm, l p = 17.5 nm, and l = 0.3 nm.The influences of the contour length of the PEG chain and the PA rods are studied.First, we investigate the changes in stiffness as we modify the ratio between the contour length of a PEG chain and that of a PA rod.Figures 9a and 9b   content, defined by the volumetric deformation J, respectively.
In the case of a constant ρ s , the volumetric ratio J depends on the contour length of the PEG chains.Maintaining a constant volumetric deformation J refers to a network in which the initial end-to-end distance of the chains depends on the number of repeat units in the PEG chain, thereby tuning the initial elasticity of the network.As expected, in the limit L PEG /L PA → ∞, corresponding to PEGDA hydrogels with PEG chains that are significantly longer than the PA rods, the limit E = E 0 is recovered.As this ratio decreases, the stiffness of the material increases due to the stiffness of the PA rods and the PA−PA interactions.
We also find that for a given ratio L PEG /L PA , increasing the contour length of the PEG chains yields stiffer PEGDA networks.This trend can be explained as follows: to maintain a constant ratio between the contour length of the PEG chains and the PA rods, an increase in n requires an increase in n PA .As a result, such networks comprise longer PA rods that increase the stiffness.Interestingly, under a constant ratio ρ 0 the volumetric deformation, and consequently the water content, networks with longer PEG chains is higher.This leads to a slight softening, but the influence of the stiff PA rods dominates.Considering a constant volumetric deformation J (and a fixed water content), longer PEG chains experience less stretch due to swelling. 37Once again, such response leads to a slight softening, but the influence of the stiff PA rods is dominant.Therefore, an overall stiffening is observed.
As previously discussed, the constitutive response to uniaxial compression is governed by the microstructural topology of the PEGDA network and the ratio L PEG /L PA .Figure 10 plots the normalized stress P x /E 0 as a function of the stretch λ for three representative ratios L PEG /L PA = 0.7, 1, and 2. As expected, decreasing the ratio L PEG /L PA results in stiffer PEGDA hydrogel networks.As previously mentioned, the response of networks characterized by L PEG /L PA = 2 is mostly governed by the entropic deformation of the PEG chains while the constitutive behavior of hydrogels with L PEG /L PA = 0.7 is governed by PA−PA interactions.

DISCUSSION AND CONCLUSIONS
In this work, we propose a microstructurally motivated model to explain the origin of the properties and response of PEGDA hydrogels.The unique structure of these gels, which comprise stiff PA rods that interconnect entropic PEG chains, leads to an interesting overall mechanical behavior that is governed by three main mechanisms: (1) the entropic response of the PEG chains, (2) the deformation of the PA rods, and (3) the PA− PA interactions.When subjected to an external loading, the PEG chains elongate or shorten entropically and apply forces on the PA rods, which in turn experience a slight bending.By considering the PEG chains as freely jointed, expressions for the energy associated with these two phenomena is derived.The PA−PA interactions are also crucial to the overall response.Specifically, the deformation of PEGDA networks with short PEG chains and long PA rods is dominated by the mechanical restrictions imposed by the rods.Such restrictions include limited rotational freedom which follows steric interactions due to the anisotropic nature of the PA rods.
The framework results in a model that can be implemented numerically to capture the mechanical response and the microstructural evolution of the network during loading.In addition, closed-form approximations for the Young's modulus and the stress that develop in networks with long chains and relatively short PA rods are developed.
To validate the model, we fabricated nine experimental hydrogel samples with different molecular weights and concentrations and subjected them to uniaxial compression.The model is found to be in agreement with the experimental findings, thus strengthening its merit.The closed-form solutions provide excellent accuracy for small to moderate deformations with a maximum error of <10% for large compressive in PEGDA networks with low molecular weight.
Our findings show that in hydrogels networks prepared with higher molecular weight PEGDA, the PEG chains are long while the PA rods are short, and accordingly the stiffness and the overall constitutive response is governed by the entropic stretching of the PEG chains.Networks with lower molecular weight PEGDA are characterized by long PA rods and relatively short PEG chains, resulting in a stiffer hydrogel with a constitutive behavior that is governed by the PA rods and the PA−PA interactions.
With the aim of providing a tool for the design of PEGDA networks, we also studied the influence of the ratio L PEG /L PA between the contour length of the PEG chain and the PA rod  on the stiffness and the overall behavior.We show that in the range of L PEG /L PA ≫ 1, the hydrogel is relatively soft with a response that is mostly entropic.Networks with PEG chains that are shorter than the PA rods are typically stiffer and less extensible due to the PA−PA interactions and the stiffness of the PA rods.
In conclusion, this work sheds light on the underlying structural mechanisms that govern the response of PEGDA hydrogels.Specifically, the findings from this work can be used to program the properties and the response of these gels by appropriate microstructural modifications, including the length of the PEG chains and the PA rods, which are controlled by the molecular weight and the concentration of PEGDA during fabrication.Control over the behavior of PEGDA is crucial for many applications, including tissue engineering, drug delivery, traction force microscopy, and microfluidics.Accordingly, these insights are expected to enhance and even optimize the performance of these materials.

Figure 2 .
Figure 2. Three-dimensional schematic diagrams illustrating the differences between (a) "classical" and (b) bottlebrush hydrogel network structures.(c) Two-dimensional close-up schematic of an individual bottlebrush motif connecting a PEG chain to each PA repeat unit.PEG chains are shown in blue, and (a) PA junction points and (b, c) rod-like PA backbones in green.For clarity, (b) shows only a fraction of PEG cross-linking chains.

Figure 3 .
Figure 3. Reference and deformed state of PEGDA network.Under compression PEG chains oriented perpendicular to the loading direction entropically extend, thereby exerting forces that rotate and (slightly) bend the PA rods toward the loading direction.

Figure 4 .
Figure 4. Representative part of the molecular structure of the PEGDA-based interconnected bottlebrush network: (a) A schematic of a PA rod pointing along the N PA ( ) direction that connects n PA = 10 randomly oriented PEG chains.(b) A representative entropic force from the PEG chain along the PA rod f N ( ) c K k PA ( , , ) First, consider the N PA (α) PA chains in the PEGDA network which are aligned along the direction N PA ( ) in the reference configuration.Here, the superscript α is used to denote a specific PA direction.It is convenient to define a local c o o r d i n a t e s y s t e m unit vectors that span

Figure 5 .
Figure 5. Illustration of a PEGDA network with (a) L PEG /L PA ≪ 1, (b) L PEG /L PA ≈ 1, and (c) L PEG /L PA ≫ 1.In the case of L PEG /L PA ≪ 1, PA rods are more likely to mechanically interact than L PEG /L PA ≫ 1.

Figure 6 .
Figure 6.Side view of a uniaxial compression test on a cylindrical sample with molecular weight M n = 2.03 kDa and concentration c = 18.4%wt/wt at the beginning (note that the uncompressed disk is 6 mm in diameter), the middle, and toward the end of the experiment.

Figure 7 .
Figure 7.Comparison of the model to the experimental findings for (a) M n = 1.07 kDa, (b) M n = 2.03 kDa, and (c) M n = 5.85 kDa.The circular markers denote the experimental findings, the continuous curves the full model, and the dashed curve the approximation in eq 27.The errors are smaller than the marker size.
plot the normalized Young's modulus E/E 0 , where E 0 = 3N PA k b T/J 1/3 is the limiting stiffness in the case of L PEG /L PA → ∞, as a function of the ratio L PEG / L PA for three representative values of n with a constant initial swollen ratio

Figure 8 .
Figure 8. Predicted Young's moduli (eq 28) as a function of the experimentally measured stiffness E for the three molecular weights M n = 1.07 kDa (triangular markers), M n = 2.03 kDa (square markers), and M n = 5.85 kDa (circular markers).The dashed line has a slope of 1 to enable comparison between the experimentally measured and the predicted stiffness.The experimental variations are smaller than the marker size.

Figure 9 .
Figure 9. Normalized Young's modulus E/E 0 as a function of the ratio L PEG /L PA for three representative values of n with (a) a constant initial swollen ratio J n / 0.3 s

Figure 10 .
Figure10.Normalized nominal stress P x /E 0 as a function of the stretch λ for three representative ratios L PEG /L PA .
Experimental video (AVI)Derivation of the energy-density of PA rods, derivation of a closed form solution for the stress associated with PA rods, experimental parameters, stress components in the PEGDA network, and comparison of the proposed model to solutions from rubber elasticity (PDF)■ AUTHOR INFORMATION Corresponding AuthorNoy Cohen − Department of Materials Science andEngineering, Technion -Israel Institute of Technology, Haifa 3200003, Israel; orcid.org/0000-0003-2224-640X;Email: noyco@technion.ac.il

Table 1 .
Summary of Model Parameters • Energy-density due to deformation of PEG chains (• = PEG), elastic bending of PA chains (• = PA), or PA−PA mechanical interactions (• = int) σ • Stress contribution due to the deformation of PEG chains (• = PEG), elastic bending of PA chains (• = PA), or PA−PA mechanical interactions (• = int) PEG /L PA → 0 J Ratio between swollen and dry volume of the PEGDA-based hydrogel network

Table 2 .
Properties of PEGDA Hydrogel Networks aQuantities marked with a # are measured directly or determined from the chemistry.M n : molecular weight of PEGDA, n: number of PEG repeat units in PEGDA; L PEG : estimated contour length of PEG chain based on a repeat unit length l ≈ 0.3 nm; c: experimental measurement of PEGDA concentration by weight; J: experimentally measured volumetric deformation; n PA : number of A repeat units in PA; L PA : contour length of a PA rod; a